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A  MOLECULAR  DYNAMICS  MODEL  OF  MELTING 
AND  GLASS  TRANSITION  IN  AN  IDEALIZED 
TWO-DIMENSIONAL  MATERIAL  -  (I) 


D.  Deng  1,  A.S.  Argon,  and  S.  Yip 
Massachusetts  Institute  of  Technology,  Cambridge,  MA  02130, 

U.S.A. 


Abstract 


^ In  a  preparatory  study  of  structural  relaxations  and  plastic  flow  in  a  two- 
dimensional  idealized  atomic  glass,  the  process  of  melting  and  quenching  through 
a  glass  transition  has  been  studied  by  computer  simulation  using  a  molecular  dy¬ 
namics  model.  In  this  model,  the  transition  from  a  solid  to  a  melt  was  observed 
to  take  place  when  liquid-like  structural  elements  composed  of  dipoles  of  5  and 
7  sided  Voronoi  polygons  percolate  through  the  two  dimensional  structure  of  dis¬ 
torted  hexagons  in  the  form  of  strings.  Such  dipoles  constitute  discrete  elements 
of  excess  free  volume  within  which  liquid  like  behaxnor  is  established  in  the  sense 
of  reduced  cohesion  or  local  elastic  moduli.  Upon  quenching  the  melt,  the  perco¬ 
lation  condition  of  liquid  like  regions  is  retained  for  under-cooled  melts  between 
the  melting  point  and  a  glass  transition  temperature  below  which  the  percolation 
condition  is  broken  and  the  thermal  expansion  is  sharply  reduced. 

c  . 

'On  leave  from  the  Institute  for  Precious  Metals  in  Kunming,  Yunnan  Province,  China. 


89  21  078 


) 


c 


r  ?  *  | 


The  simulation  which  has  used  empirical  pair  potentials  characteristic  of  &u 
and  Zt>  has  substantially  underpredicted  the  melting  and  glass  transition  tem¬ 
peratures  and  overpredicted  the  thermal  expansion  of  CukfZr^t  type  glasses. 
These  defects  of  the  model  can  be  attributed  to  the  two-dimensional  nature  of  the 
material,  which  stores  larger  concentrations  of  free  volume  than  a  correspond¬ 
ing  three-dimensional  material.  In  spite  of  these  quantitative  shortcomings,  the 
model  gives  valuable  insight  into  the  topological  features  of  the  local  atomic  con¬ 
figurations  at  melting  and  upon  vitrification. C  f  Uj\ 
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I.  INTRODUCTION 


The  modes  of  inelastic  deformation  in  amorphous  solids  on  the  atomic  scale 
below  their  glass  transition  temperatures,  i.e.,  in  their  glassy  behavior  range,  and 
how  these  are  related  to  the  structure  of  these  solids,  has  been  of  fundamen¬ 
tal  interest.  The  approaches  taken  to  elucidate  these  phenomena  in  amorphous 
metals,  polymers,  and  network  glasses  have  differed  considerably  along  the  tradi¬ 
tionally  different  views  taken  by  investigators  in  these  fields.  Nevertheless,  these 
approaches  have  in  general  divided  into  two  operationally  different  categories. 
In  the  first  category  are  those  studies  that  are  based  on  kinematically  accept¬ 
able,  strain  producing  ad-hoc  local  mechanisms  which  involve  cooperative  atom 
or  molecular  segment  motions,  and  which  are  constructed  principally  to  model 
the  inelastic  shear  resistance  [1-9].  In  the  second  category  are  mostly  studies 
based  on  detailed  computer  simulations  of  the  structure  of  atomic  [10,11]  (for 
a  review  see  [12])  and  polymeric  [13,14]  glasses,  including  some  further  analysis 
to  explore  the  mechanisms  of  inelastic  deformation  in  these  simulated  systems 
[15-17].  Related  to  the  latter  are  studies  [8,18,19]  using  analog  models  of  two- 
dimensional  atomic  glasses  based  on  the  classical  Bragg  bubble  model  [20].  The 
various  ad-hoc  models  and  computer  simulations  have  been  necessary  because 
of  the  difficulty  of  obtaining  information  by  direct  experimental  methods  on 
the  very  local  structural  alterations  that  occur  in  the  inelastic  deformation  in 
glasses.  Some  investigators  have  pursued  a  different  line  of  approach  relating 
inelastic  deformations  in  glasses  to  extensions  of  well  defined  mobile  crystal  de¬ 
fects,  such  as  vacancies  and  dislocations  (for  a  review  see  [2]).  In  our  view. 
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however,  these  have  not  furnished  many  quantitatively  useful  explanations  of 
either  the  topological  features  of  the  deformation  process  or  its  kinetics. 


In  this  and  the  accompanying  three  communications  (referred  to  as  II-IV  in 
what  follows),  we  report  the  results  of  a  series  of  computer  simulations  of  the 
structure,  its  phase  transitions,  structural  relaxations,  and  plastic  deformation 
under  stress  of  a  two-dimensional  model  atomic  glass.  In  these  simulations, 
our  goal  was  not  to  obtain  quantitatively  exact  descriptions  of  structures,  of 
physical  processes  of  phase  change,  or  of  mechanical  properties,  results  which 
can  be  compared  directly  with,  experimental  data  on  metallic  glasses.  There 
are  certain  topological  limitations  of  a  two-dimensional  simulation  which  does 
not  permit  this.  Rather,  our  goal  has  been  to  elucidate  complex  topological 
processes  of  a  cooperative  nature,  and  to  develop  qualitative  understanding  and 
semi-quantitative  scaling  laws  for  them. 

One  principal  objective  of  the  study  was  the  simulation  of  structural  aging 
and  the  subsequent  cooperative  atom  motions  occurring  in  local  shear  transfor¬ 
mations  upon  large  strain  shearing.  To  obtain  the  initial  structures  for  these 
studies  in  glassy  assemblies  in  a  non-arbitrary  way,  it  was  decided  to  obtain  all 
such  initial  states  by  melting  and  quenching  the  two-dimensional  atom  assem¬ 
blies.  In  this  communication,  we  discuss  only  our  findings  related  to  melting 
and  the  glass  transition  upon  quenching. 

In  several  instances,  we  have  encountered  significant  constraints  on  the  kine¬ 
matics  and  the  distributed  nature  of  the  processes  of  shear  relaxation  even  in 
relatively  large  two-dimensional  periodic  cells,  indicating  the  wisdom  of  pre¬ 
ferring  a  large  two-dimensional  array  of  a  fixed  number  of  atoms,  rather  than 
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distributing  them  in  a  more  confined  but  three-dimensional  periodic  cell  -  when 
the  overall  number  could  not  be  arbitrarily  increased  for  reasons  of  tractability 
in  the  simulation  -  not  to  mention  the  difficulty  in  visualization  of  the  local 
transformations. 

II.  DETAILS  OF  THE  SIMULATION 

2.1  The  Atomic  Potentials 


The  simulations  reported  here  and  in  the  accompanying  communications 
make  use  of  standard  molecular  dynamics  (MD)  simulation  techniques  [21].  In 
(MD)  Newton’s  equations  of  motion  for  a  two-dimensional  system  of  N,  par¬ 
ticles  are  integrated  using  a  fifth-order  predictor-corrector  algorithm.  Periodic 
boundary  conditions  are  imposed  on  the  simulation  cell,  and  through  a 
Lagrangian  formulation  [22],  the  components  of  the  basis  vectors  defining  the 
cell  are  allowed  to  vary  in  accordance  to  any  imbalance  between  the  internal 
stress  and  an  externally  prescribed  stress.  The  system  temperature  is  defined 
in  terms  of  the  average  kinetic  energy  of  a  particle.  To  maintain  the  system 
at  a  particular  temperature,  the  particle  velocities  are  rescaled  at  every  time 
step  of  the  simulation.  If  the  system  has  reached  thermal  equilibrium,  then 
the  instantaneous  temperature,  without  velocity  rescaling,  will  fluctuate  about 
a  mean  value  which  does  not  change  in  time. 
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The  atomic  interactions  in  the  present  simulations  are  given  by  pair  poten¬ 


tials  of  the  Lennard- Jones  type  1  [10]: 


*(r«)  =  E, 


-  aVjj  +  b' 


(1) 


where  $(rtJ)  is  the  potential  energy  of  interaction  between  atoms  i  and  j  at 
a  separation  distance  r„,  r0  and  E0  are  length  and  energy  parameters  of  the 
potential,  and  the  constants  a’  and  b'  are  chosen  to  obtain  smooth  truncation  of 
the  potential  some  place  between  the  third  and  fourth  nearest  neighbors.  The 
exponents  m  and  n  are  phenomenologically  assigned  to  obtain  a  desired  form 
for  the  potential. 

In  our  simulations,  two  different  forms  of  the  potential  were  used.  In  the 
early  phase  of  the  study  of  melting  and  stability  of  the  quenched  glassy  state  of  a 
single  component  solid,  a  Lennard-Jones  potential  was  used  with  m  =  12,  n  =  6, 
r0  =  3.405.A,  and  an  energy  to  atomic  mass  ratio  E0/ma  =  2.485xl08cm2/sec2, 
to  give  a  fundamental  atomic  period  rj  =  r0\Jm0/E0  =  2.16xl0~12sec 
(E0  =  1.034xl0-2eV,  and  m0  =  6.656xl0-23g).  These  particular  parameter 
values  are  appropriate  to  the  noble  gas  solid  argon.  The  simulations  were  carried 
out  under  two  constant  external  pressures:  p*  =  0,  and  1,  in  units  of  ir  =  E0/r\ 
(=  41.9  MPA  in  3-D),  to  explore  the  effect  of  pressure  on  structural  stability. 
In  a  typical  simulation  run,  the  time  step  size  has  been  taken  at  1%  of  the 

fundamental  period  rj,  which  makes  it  2.16xl0-usec  in  this  case.  As  we  will 
1For  nomenclature  of  symbols,  refer  to  the  table  at  the  end  of  this  paper. 


4 


discuss  below  briefly,  this  single  component  solid  crystallized  exceedingly  rapidly. 
Therefore,  simulations  with  this  system  were  not  carried  much  further. 

In  the  main  simulation,  which  will  be  the  subject  of  this  and  the  accompany¬ 
ing  communications,  a  modified  Lennard-Jones  potential  of  the  type  introduced 
by  Kobayashi  et  al.  [10]  for  a  two-component  glass,  such  as  CuxZrx.x,  was  used. 
For  the  system  considered  here,  x  =  0.5,  m  =  8,  n  =  4.  The  Cu  and  Zr  atoms 
have  been  labelled  as  A  and  B  respectively,  and  their  characteristic  parameters 
have  been  taken  as  follows:  r0  —  2.556.A,  m0  =  10.63xl0~”g,  E0  =  0.15eV  for 
the  Cu  (A)  atoms,  giving  a  fundamental  atomic  period  rjA  =  5.4xl0~issec,  a 
fundamental  stress  unit  irA  =  1.437GPa  (in  3-D).  The  interactions  between  the 
Zr  atoms  (B-B)  and  the  Cu  and  Zr  atoms  (A-B)  were  obtained  according  to 
the  following  ratios:  roBBlr0  =  1.244;  roABjr0  =  0.5  (1  +  roBB/r0)  =  1.122; 
EoBB/E0  =  1.94;  EoAB/E0  =  E0bb/E0)1/ 2  =  1.393.  Furthermore,  the  two 
coefficients  a'  and  b'  were  chosen  such  that  interatomic  forces  vanish  at  a  critical 
cut-off  distance  of  rc  =  2.5  r0.  The  resulting  pair  potential  for  copper  is  shown 
in  Fig.  1.  The  two  other  pair  potentials  look  very  similar.  To  enhance  stability 
in  this  simulation,  an  external  pressure  of  p/?r  =  1.0  was  chosen  (i.e.,  1.437  GPa 
in  3-D).  This  is  a  substantial  pressure,  and  its  effect  on  the  quantitative  aspects 
of  the  simulation  must  be  kept  in  mind.  The  time  step  size  was  taken  at  5%  of 
the  fundamental  period,  i.e.,  at  l.lxlO-14  sec. 

The  method  by  which  the  simulation  cell  is  allowed  to  change  volume  and 
shape  is  based  on  a  Lagrangian  which  contains  a  kinetic  energy  term  associated 
with  the  motions  of  the  cell  vectors  [22,23].  A  mass  parameter  appears  in  this 
term;  whiie  it  can  be  assigned  any  value,  a  proper  choice  is  important  because 
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it  governs  how  quickly  the  cell  vectors  respond  to  the  imbalances  between  the 
prescribed  external  pressure  and  the  fluctuating  internal  pressure.  We  discuss 
the  effect  of  this  so-called  wall  mass  on  the  simulations  and  its  proper  choice  in 
Appendix  I  of  the  accompanying  paper  on  structural  relaxations,  which  we  will 
refer  to  here  as  III  [30].  The  magnitude  of  this  mw  in  the  present  simulation 
was  4. 

2.2  Simulation  Cell 

The  simulation  cells  that  were  chosen  in  this  study  contained  144  atoms 
initially  arranged  on  a  12x12  grid  in  a  hexagonal  close  packed  pattern  (see  Fig. 
9a).  Because  of  the  hexagonal  symmetry,  the  simulation  cell  was  11.5%  longer 
in  the  horizontal  direction  than  in  the  vertical  direction.  Areas  allocated  to 
each  atom  were  delineated  by  Voronoi  polygons  constructed  in  the  conventional 
way.  This  was  done  under  every  circumstance,  and  proved  to  be  of  great  value 
in  visualization  of  distortions,  but  also  in  determination  of  strain  increments  at 
each  atomic  site  and  the  level  of  coordination  of  atoms  to  their  surroundings.  To 
simulate  properties  of  a  two-dimensional  solid  of  infinite  extent  in  the  plane  of 
simulation,  periodic  boundary  conditions  were  used  throughout  the  simulations. 
In  the  simulations  of  the  binary  atom  solids  of  Cu0.sZro.s ,  the  Cu(A)  and  Zr( B) 
atoms  were  initially  assigned  randomly. 
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S.S  Melting  and  Quenching  Simulations 

The  melting  of  the  two-dimensional  crystal  lattice  was  achieved  by  step-wise 
increase  of  the  system  temperature.  The  dimensionless  temperature  T*{=  T/6 ) 
was  increased  in  steps  that  ranged  from  0.025  to  0.1  (the  normalized  temperature 
6  =  E0/k  is  120K  for  the  one-component  solid  and  1739K  for  the  two-component 
Cu0iZr0i  solid).  Between  increments,  the  temperature  was  held  constant  by 
rescaling  the  particle  velocities.  To  achieve  equilibrium  as  defined  by  cumulative- 
averaged  properties  no  longer  changing  with  time,  it  was  necessary  to  equilibrate 
for  3  —  4x10s  time  steps  at  temperatures  below  the  melting  temperature  and  for 
3  — 5x10s  time  steps  above  the  melting  temperature.  During  the  melting  process, 
where  large  and  discontinuous  changes  occur  in  the  system,  often  more  time  steps 
were  necessary  to  achieve  stabilized  melt  structures.  All  temperature  dependent 
intensive  properties  of  the  equilibrated  system  for  the  heating  runs  have  been 
calculated  after  structural  stability  was  attained  at  the  new  temperature,  judged 
by  an  absence  of  further  change  in  the  average  volume  per  atom.  The  same 
was  done  on  cooling  in  the  melt,  down  to  where  time  dependent  structural 
relaxations  take  progressively  longer  periods  of  time.  In  this  case,  the  structural 
relaxations  were  only  carried  out  partially  for  each  cooling  step  to  freeze  in  a 
structure  typical  of  that  of  the  subcooled  melt.  In  this  range,  each  temperature 
reduction  of  about  0.05  -  0.1  was  followed  by  only  10s  time  steps  of  stabilization. 
Such  quenching  simulations  amount  to  very  high  rates  of  cooling  in  the  range 
of  3xlOnK/sec  for  the  one-component  system  and  3xl010K/sec  for  the  two- 
component  systems. 
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The  atom  coordinates  for  each  state  in  the  melting  and  quenching  runs  were 


stored  for  future  use  in  relaxation  studies. 

III.  RESULTS 


S.l  Temperature  Dependent  Changes  in  Atomic  Volume,  Enthalpy, 
and  Potential  Energy 


The  changes  in  the  average  atomic  volume  (average  Voronoi  area)  Cl*,  en¬ 
thalpy  per  atom  H* ,  and  potential  energy  per  atom  V’,  all  in  normalized  units, 
resulting  from  the  changes  in  system  temperature  have  been  determined  for  the 
melting  and  quenching  runs  for  both  the  single-component  solid  under  two  pres¬ 
sures  (0,  1.0),  and  the  two-component  solid  under  a  pressure  of  1.0,  for  each 
step  of  the  runs. 

The  problem  of  relating  the  two-dimensional  (2-D)  results  to  three- 
dimensional  (3-D)  properties  has  been  dealt  with  systematically  as  follows: 

a.  Stresses  and  pressures  are  defined  for  use  in  (3-D),  but  operationally  are 
handled  in  intermediate  level  computations  in  (2-D),  i.e.,  as  forces  per 
peripheral  length.  In  final  presentations,  all  results  are  normalized  by  the 
three-dimensional  pressure  unit  ( E0/r\ )  so  the  reader  can  determine  (3-D) 
properties  if  desired. 
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b.  All  intensive  properties  of  interest  obtained  from  the  simulation  have  been 
presented  per  atom  when  given  in  dimensioned  units.  In  all  other  cases, 


they  have  been  given  in  non-dimensional  form,  but  still  per  atom. 

Thus,  the  (2-D)  simulation  cell  has  been  treated  for  (3-D)  purposes,  as  if  it 
had  a  thickness  of  r0,  although  in  every  instance,  the  simulation  was  constrained 
to  lie  in  a  plane. 

While  the  average  Voronoi  polygon  area  per  atom  fl*  in  units  of  the  funda¬ 
mental  area  r \  is  obtained  directly  from  the  simulation  plane,  the  instantaneous 
potential  energy  per  atom  V*  and  enthalpy  per  atom  ff',  on  a  per  atom  basis 
were  obtained  as  follows: 


V 


1 

nE0 


M. 


H* 


1 

nE0 


rriivl  +  V*  +  pfi* 


(3) 


The  results  for  H*,  H *,  and  V *  are  plotted  as  a  function  of  T*  in  Figs.  2 
and  3  at  two  external  pressures  p’  of  0  and  1.0  (in  units  of  n  =  E0/rl)  for  the 
single-component  solid,  and  Fig.  4  for  the  two-component  solid  (constructed  of 
equal  numbers  of  atoms  A  and  B)  at  an  external  pressure  of  1.0. 
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In  ail  cases,  the  properties,  ft*,  H *,  and  V*  increase  monotonicaily  with 
increasing  temperature.  The  changes  in  the  properties  of  the  one-component 
solid  (Cases  I  and  II),  and  the  two-component  solid  (Case  III)  with  increasing 
temperature  are  summarized  in  Table  I,  where  a0  and  am  are  the  coefficients 
of  linear  thermal  expansion  at  zero  temperature  and  the  melting  point  respec¬ 
tively.  It  can  be  seen  that  the  ratio  am/ a0  is  6-9  for  the  one  component  solid  and 
2.5  for  the  two-component  solid.  Across  the  melting  point,  a  does  not  increase 
significantly.  The  same  is  true  for  the  rate  of  change  of  enthalpy  and  poten¬ 
tial  energy  per  atom.  This  suggests  that  the  melting  observed  here  is  more  a 
process  of  topological  transition,  and  manifestations  of  a  discontinuous  behavior 
in  the  rates  of  change  of  intensive  properties  are  weak.  This  conclusion  will 
be  reinforced  later  by  inspection  of  the  topological  features  associated  with  the 
melting  process.  The  latter  are  to  be  expected  because  the  limited  system  size 
and  simulation  durations  make  the  observed  transition  necessarily  smeared  out. 

The  fractional  linear  expansion  in  these  model  solids  upon  melting  range 
from  about  2.14xl0~s  for  the  one-component  system  to  5.49xl0_s  for  the 
two-component  system,  giving  projected  volumetric  expansions  on  melting  of 
6.42xl0"1  for  the  one-component  solid  and  1.65xl0_l  for  the  two-component 
solid.  These  results  do  not  compare  too  well  with  the  known  volume  expansions 
of  noble  gas  solids  (12.7%  for  Ar)  and  close-packed  metals  (4.2%  for  Cu). 

The  absolute  agreement  in  the  melting  points  is  also  unsatisfactory.  In  the 
one-component  solid,  the  melting  point  under  no  pressure  is  42K  and  under  a 
pressure  of  p*  =  1.0  (42  MPa)  is  72K.  Compared  with  the  melting  point  of  solid 
Argon  of  84K  at  atmospheric  pressure,  both  values  are  too  low.  The  same  is  true 
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for  the  two-component  solid,  which  at  a  pressure  of  p*  =  1.0  (1.437  GPa),  has 
a  melting  point  of  434K,  while  the  melting  point  of  the  intermetallic  compound 
Cxio.sZro.h  which  it  is  to  simulate,  has  a  melting  point  of  1253K  at  atmospheric 
pressure. 

It  is  worthwhile  to  note  here  that  the  phase  diagram  of  the  two-dimensional 
Lennard-Jones  system  has  been  studied  by  Monte  Carlo  simulation  [24].  The 
triple  point  conditions  are  known  to  be  p\  —  .815,  =  .415,  and  p*t  =  .0056. 

Since  p“t  is  essentially  zero,  we  can  compare  the  melting  temperature  T £  and 
corresponding  density  p'm  in  Fig.  2  with  these  values.  We  find  that  =  .353 
and  p*m  —  .812.  We  believe  that  the  reason  why  our  value  is  lower  than 
the  literature  value  is  because  our  range  cutoff  re  is  somewhat  shorter,  and  this 
effect  is  known  to  lower  the  melting  point.  In  addition,  we  have  not  taken  into 
account  the  long  range  correction  to  the  pressure  or  correction  for  system  size 
dependence. 

The  pressure  dependence  of  the  melting  point,  as  simulated  in  our  compu¬ 
tations  for  the  one-component  solid  is  11.3K/GPa.  This  is  difficult  to  compare 
with  actual  experimental  results,  since  the  corresponding  information  is  un¬ 
available  for  noble  gas  solids.  On  a  normalized  basis,  however,  a  comparison 
is  possible  where  the  result  from  the  simulation  ( ATm/Tm)/(AP/n )  is  0.526. 
The  corresponding  information  for  an  alkali  metal,  such  as  Na  is  0.338,  where 
7r  =  ( E0/rl )  is  taken  as  0.65  GPa,  as  scaled  down  from  the  given  ratio  of  Cu,  on 
the  basis  of  their  different  elastic  properties.  Clearly,  the  corresponding  effect 
for  a  noble  gas  solid,  such  as  Argon,  should  be  larger  and  more  comparable  with 
our  simulation. 
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From  the  information  discussed  above,  we  conclude  that  our  simulation  in 
2-D  underestimates  the  absolute  melting  points,  but  gives  qualitatively  accept¬ 
able  results  for  the  dependence  of  the  melting  point  on  pressure  on  a  normalized 
basis.  Furthermore,  the  results  of  the  simulation  on  the  linear  coefficient  of 
expansion  for  the  two-component  solid  at  7.35xlO~5K-1  at  low  temperatures  is 
quite  high  when  compared  with  what  might  be  expected  for  Cuo.sZr0.s,  esti¬ 
mated  to  be  about  l.lxlO-5K_1  as  an  average  for  Cu  and  Zr ,  at  equal  volume 
fractions.  This  behavior  of  the  simulation,  however,  is  not  unexpected,  since 
2-D  models  of  disordered  media  are  expected  to  give  overestimates  of  specific 
volume  and  related  properties,  such  as  thermal  expansion.  Since  the  purpose 
of  the  simulation  of  melting  was  principally  to  obtain  non-arbitrary  starting  co¬ 
ordinates  for  atoms  in  a  disordered  medium  in  the  solid  state  upon  quenching 
from  the  molten  state,  we  do  not  attach  too  much  importance  to  the  relatively 
poor  agreement  between  the  results  of  the  simulation  and  the  absolute  melting 
behavior.  We  will,  however,  find  many  points  of  interest  in  the  topological  fea¬ 
tures  of  the  melt  and  its  alterations  upon  quenching,  as  we  will  discuss  these 
results  below. 

When  the  melt  is  cooled  by  systematically  reducing  the  temperature  in  steps 
with  intervening  short  periods  of  structural  relaxation  to  achieve  a  condition  of 
very  slow  change  of  properties  (c.a.  1000  time  steps),  the  intensive  properties 
decrease  continuously  and  monotonically.  While  the  decreases  above  the  melting 
temperature  follow  the  same  behavior  of  the  heating  curve,  in  reverse,  this  is 
not  so  below  the  melting  point,  as  Figs.  2-4  clearly  show.  At  these  cooling  rates, 
which  are  in  the  range  of  10nK$ec-1,  the  liquid  phase  has  no  time  to  crystallize 
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by  nucleation,  but  follows  a  temperature  dependent  contraction  roughly  parallel 
to  that  of  the  solid  hexagonal  close-packed  phase,  in  an  undercooled  state,  for 
the  one  and  two  component  materials  under  pressure  (Figs.  3  and  4).  In  the 
one-component  material,  the  slope  of  the  cooling  curve  always  remains  larger 
than  that  of  the  corresponding  crystal,  making  the  two  curves  meet  before  reach¬ 
ing  zero  temperature.  This  indicates  that  the  one-component  system  under  no 
pressure  in  the  undercooled  state  is  unstable  and  crystallizes  even  at  these  very 
high  rates  of  quenching.  In  comparison,  both  the  one-component  system  under 
pressure,  as  well  as  the  two-component  system  under  pressure,  do  not  join  onto 
the  hexagonal  crystal  and  remain  at  least  quasi-stable  in  the  disordered  form. 

In  real  systems,  when  a  changeover  occurs  in  the  rates  of  change  of  intensive 
properties  in  a  very  narrow  range  of  temperature  from  a  liquid-like  behavior  to 
solid-like  behavior,  without  any  discontinuity  in  the  property  itself,  the  material 
is  said  to  undergo  a  glass  transition.  This  is  indicated  fairly  clearly  in  Fig.  4  for 
the  two-component  material,  where  a  glass  transition  is  found  at  a  normalized 
temperature  of  0.18,  in  comparison  to  the  melting  point  at  0.25. 

There  exist  several  simulation  studies  of  the  Lennard-Jones  system  adopted 
here.  Monte  Carlo  results  have  been  reported  on  melting  in  two  dimensions 
[25],  and  on  the  glass  transition  in  two  dimensions  in  one  and  two  component 
systems  [26].  In  three  dimensions,  molecular  dynamics  results  on  melting  [27], 
crystallization  [27],  and  glass  transition  [28]  are  available. 


S.2  Changes  in  Topological  Features 


3.2.1  Radial  distribution  function  (RDF) 


The  changes  in  the  RDF  of  the  three  materials  that  occur  upon  heating 
are  quite  similar.  In  Figs.  5a-5d,  we  show  the  typical  sequence  for  the  two- 
component  solid.  At  increasing  temperatures,  the  peaks  occur  at  larger  distances 
and  they  become  broader  and  less  intense. 

Two  important  qualitative  changes  occur  in  the  RDF  near  the  melting  point. 
First,  the  second  and  third  peaks,  which  are  clearly  separated  in  the  expanded 
hexagonal  close-packed  structure  below  Tm  merge  together  and  form  a  very 
broad  and  short  peak  at  T  >  Tm.  This  also  occurs  for  the  fourth  and  fifth 
peaks,  which  also  spread  out  and  fuse  together  in  the  melting  process.  Second, 
the  clear  separation  between  peaks  present  below  Tm  begins  to  be  filled  up  near 
the  melting  point.  This  is  particularly  evident  in  the  larger  separation  between 
the  first  and  second,  and  between  the  third  and  fourth  peaks.  Thus,  above  the 
melting  point,  the  RDF  has  non-zero  values  at  all  separations  and  the  previously 
well  defined  first  six  peaks  fuse  into  three  broad  and  spread-out  peaks.  Some  of 
these  changes  are  also  observable  just  below  the  melting  point.  This,  as  we  have 
already  pointed  out  is  a  result  of  the  limited  system  size  and  short  durations  of 
the  simulation  which  spreads  out  transitions. 

The  changes  in  the  RDF  at  melting  reflect  the  changes  and  smearing  out  of 
structural  order.  This  is  best  followed  from  the  ideal  (2-D)  hexagonal  lattice 
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shown  in  Fig.  6  and  its  RDF  reproduced  immediately  under  it.  Clearly,  the 
second  and  third,  and  the  fourth  and  fifth  peaks  are  relatively  close  together, 
so  that  increased  amplitude  of  thermal  motion  in  the  solid  stat.  anu  ventually 
the  even  larger  amplitudes  of  atom  motions  possible  in  the  molten  state  should 
smear  out  and  fuse  these  peaks  together. 

The  changes  in  the  RDF  upon  quenching  for  the  three  materials  are  quali¬ 
tatively  different  in  the  one-component  solid  from  the  two-component  solid  in 
the  time  range  of  the  simulation.  The  one-component  solid  under  no  confining 
pressure  crystallizes  nearly  completely,  as  Figs.  7a-7d  show  clearly.  The  same 
is  also  true  under  even  a  pressure  of  1.0.  In  comparison,  for  the  same  quench¬ 
ing  schedule,  the  two-component  material  under  a  pressure  of  1.0  is  much  more 
stable  and  reaches  a  reasonably  stable  state  at  a  temperature  of  0.06,  as  Figs. 
8a-8c  show.  Even  here,  however,  as  is  seen  in  Fig.  8c,  the  second  and  third 
peaks  have  clearly  begun  to  separate,  and  there  are  indications  that  the  fifth 
peak  is  about  to  separate  from  the  fourth  peak.  Nevertheless,  the  distribution 
is  still  continuous  and  without  any  breaks.  We  take  this  instantaneous  state 
of  well  relaxed  structure  as  a  borderline  case  between  a  defective  crystal  and  a 
glass  with  considerable  topological  short  range  order. 

3.2.2  Topological  changes  in  melting 


The  topological  changes  occurring  in  melting  for  the  one-component  and  two- 
component  systems  under  a  normalized  pressure  of  unity  (p*  =  1)  are  shown  in 
Figs.  9  and  10.  The  areas  per  atom  are  delineated  by  Voronoi  polygons.  Their 
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changing  sizes,  shapes,  and  associations  reveal  important  details  which  cannot 
be  perceived  from  the  RDF. 

In  the  sequences  in  Fig.  9,  the  changes  in  the  one-component  material  are 
shown  with  increasing  temperature.  It  must  be  born  in  mind  that  while  each 
picture  reflects  the  structure  in  near  equilibrium  after  1  —  4x10s  time  steps,  they 
are  nevertheless  only  a  snap  shot  in  time,  showing  only  one  of  a  large  number  of 
statistically  equivalent  configurations  with  the  same  overall  energy.  Figure  9a 
shows  the  simulation  cell  of  1 44  atoms,  which  still  shows  a  relatively  well-ordered 
hexagonal  atomic  arrangement  of  the  solid  at  half  its  melting  point  (T*  =  0.3). 
While  all  polygons  are  still  hexagons  at  this  temperature,  careful  inspection  of 
these  shows  that  some  of  them  have  acquired  thermal  motion  induced  small 
distortions.  At  a  normalized  temperature  of  0.5  ( T  =  0.83Tm),  all  polygons  are 
still  hexagons,  but  many  of  them  show  quite  large  distortions.  In  Fig.  9b,  at  the 
melting  point  (T„  =  0.6)  after  only  1000  time  steps  before  complete  melting  has 
occurred,  qualitatively  different  topological  changes  are  apparent.  Now,  most  of 
the  hexagons  are  quite  visibly  distorted,  and  many  of  the  hexagons  have  been 
transformed  into  pentagons  (5)  and  heptagons  (7).  Furthermore,  accompanying 
the  5  and  7  sided  polygons,  corners  sharing  4  edges  also  appear  as  a  transient 
form  for  the  first  time.  What  is  more,  the  5  and  7  sided  polygons  show  an  affinity 
to  each  other.  Figure  9c  shows  the  structure  at  the  melting  point  after  4000  time 
steps  where  structural  equilibrium  has  been  reached  and  no  further  expansion 
occurs.  The  simulation  cell  in  the  molten  state  has  acquired  a  significant  shear 
distortion  without  any  applied  shear  stress.  This  cam  be  traced  to  an  artifact 
related  to  the  smallness  of  the  simulation  cell.  In  a  cell  of  amorphous  material 
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of  only  144  atoms  subjected  to  a  periodic  boundary  condition,  the  repetition  in 
space  of  the  borders  of  the  cell  constitutes  a  paracrystal  in  which  the  periodic 
border  material  makes  up  a  significant  volume  fraction.  This  volume  fraction 
can  increase  even  further  if  the  cell  is  allowed  to  shear  in  the  liquid  state,  where 
the  internal  resistances  to  deformation  have  become  very  small.  We  ignore  this 
distortion.  More  importantly,  in  the  simulation  cell  of  the  molten  material,  a 
large  number  of  the  Voronoi  polygons  have  now  become  5  and  7  sided,  indi¬ 
cated  by  the  numbers  in  the  polygons,  and  these  have  associated  themselves  in 
strings,  in  which  4  edged  corners  also  appear  frequently.  In  most  instances,  the 
strings  percolate  through  the  cell,  as  is  shown  in  the  bottom  of  Fig.  9c.  At  still 
higher  temperatures,  the  percolating  string-like  association  of  the  5  and  .7  sided 
polygons  in  the  simulation  cell  become  an  even  more  prominent  feature. 

The  corresponding  distortions  upon  melting  in  the  two-component  material 
under  a  normalized  pressure  of  unity  ( p *  =  l)  are  shown  in  the  sequences  in  Fig. 
10.  In  Fig.  10a,  the  initial  hexagonal  arrangement  is  shown  at 
T‘  =  0.1  [T  =  0.4Tm).  The  Voronoi  polygons  containing  the  stiffer  Zr  atoms 
are  identified  by  the  small  circles  in  the  polygons.  Because  of  the  differences  in 
size  and  local  binding  between  the  Cu  and  Zr  atoms,  many  polygons  appear 
already  significantly  distorted.  Any  apparent  clustering  of  the  Zr  atoms  in  the 
simulation  cell  is  illusory,  since  the  Cu  and  Zr  atoms  were  assigned  to  polygons 
purely  randomly.  At  T*  =  0.2(=  0.8Tm),  again  the  hexagons  show  similarly 
increased  distortions,  as  in  the  one-component  material.  Figure  10b  shows  the 
high  temperature  melt  at  T *  =  0.3 (T  =  1.2 Tm)  with  percolating  strings  of  5 
and  7  sided  polygons.  At  the  same  homologous  temperature  with  respect  to  the 
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melting  point,  the  departures  from  regularity  of  the  hexagons  are  larger  in  the 
two-component  system  than  in  the  one-component  system. 

We  will  demonstrate  below  that  the  5  and  7  sided  polygons  usually  associate 
themselves  into  a  structural  dipole  of  liquid-like  material,  and  that  the  further 
percolating  string-like  associations  of  these  dipoles  imparts  fluid-like  behavior 
to  the  structure. 

3.2.3  Topological  changes  in  undercooled  melts 
and  in  a  glass  transition 

The  topological  changes  in  the  arrangement  of  atoms  undergoing  quenching 
is  shown  in  Fig.  11  for  the  two-component  system,  which  is  typical  for  all 
systems  that  were  studied,  including  the  one-component  system.  Figure  11a 
shows  the  two-component  system  at  T*  =  0.3(T  =  1.2 Tm)  in  the  melt.  A  large 
and  well-formed  percolation  string  of  5  and  7  sided  polygons  is  clearly  visible. 
There  are  other  broken  strings  of  liquid-like  polygon  dipoles  also  visible,  which 
can  connect  up  to  the  main  chain  and  thereby  replace  other  segments  of  that 
chain  by  a  statistically  equivalent  chain  in  other  snapshots  in  time.  Figure  lib 
shows  the  polygons  in  the  simulation  cell  at  T*  =  0.2(7’  =  0.8Tm),  which  is  an 
undercooled  melt  between  Tm  and  Tt.  A  percolating  string  of  liquid- like  material 
is  still  evident  in  the  figure.  It  is  well  to  recall  that  the  simulation  involving  step¬ 
wise  lowering  of  the  temperature  followed  by  only  10s  time  steps  for  equilibration 
of  the  structure  constitutes  a  quenching  rate  of  about  10n/fsec_1.  In  Fig.  11c, 
the  system  is  shown  at  T *  =  0.14(=  0.78T,  =  0.56 Tm),  where  many  associated 


dipoles  of  5  and  7  sided  polygons  are  visible,  but  the  string-like  percolation 
of  these  through  the  simulation  cell  is  now  broken.  Thus,  the  material  has 
apparently  regained  its  solid-like  behavior.  Finally,  in  Fig.  lid,  at  a  temperature 
of  T'  =  0.06(=  0.33T,  =  0.24Tm),  some  5  and  7  sided  polygon  dipoles  still  exist, 
but  are  now  isolated  in  a  background  of  short  range  ordered  material.  This 
represents  the  well-relaxed  glassy  state  illustrated  by  the  RDF  in  Fig.  8d. 


IV.  DISCUSSION 


4.1  Melting  and  Glass  Transition 

Our  simulation  of  two-dimensional  hexagonal  crystals  has  shown  that  in  the 
first-order  transition  of  melting,  the  system  does  not  expand  by  spreading  the 
disorder  quasi-uniformly,  but  rather  acquires  well-defined  imperfections  having 
excess  free  volume.  The  principal  ingredient  of  these  imperfections  are  structural 
dipoles  of  associated  Voronoi  polygons  of  5  and  7  sides. 

When  the  system  is  molten  so  that  solid-like  behavior  is  replaced  by  liquid¬ 
like  behavior,  a  critical  concentration  of  the  material  must  acquire  a  radically 
lowered  level  of  cohesion.  This  occurs  through  the  establishment  of  a  percola¬ 
tion  condition,  in  which  strings  of  5  and  7  sided  polygon  dipoles  pervade  the 
entire  simulation  cell  in  one  or  more  directions  and  hence,  through  the  periodic 
boundary  conditions,  reach  over  the  entire  infinite  (2-D)  space. 

Clearly,  for  this  picture  to  have  meaning,  the  5-7  dipole  must  have  elemental 
liquid-like  properties.  That  this  is  so,  is  shown  in  Fig.  12,  which  gives  the  average 
areas  of  the  Voronoi  polygons  in  units  of  r\,  of  5, 6,  and  7  sides,  as  measured  from 
actual  simulation  cells  at  a  temperature  of  T*  =0.1.  These  relative  differences 
do  not  change  much  with  temperature  up  to  the  melting  point.  As  the  figure 
shows,  the  average  heptagon  is  about  10%  larger  than  the  average  hexagon, 
which  on  the  average  is  only  3%  larger  than  the  average  pentagon.  Thus,  the 
area  of  the  average  dipole  of  5-7  sided  Voronoi  polygons  is  about  7%  larger  than 
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a  pair  of  average  hexagons.  Since  the  average  hexagons  themselves  are  about 
9%  larger  at  the  melting  point  than  the  hexagons  at  0 °K,  the  average  dipole  in 
the  melt  is  about  16%  larger  than  a  pair  of  hexagons  in  a  (2-D)  material  at  rest 
at  OK.  This  is  a  substantial  increase,  which  should  put  the  material  of  the  dipole 
at  or  beyond  the  decohesion  strain,  and  give  it  indeed  liquid-like  properties.  We 
show  in  the  accompanying  communications  [29,31]  (to  be  referred  to  here  as  (II) 
and  (IV))  that  both  the  structural  relaxations  and  plastic  deformation  at  low 
temperatures  are  strongly  influenced  by  the  liquid-like  material. 

The  importance  of  such  liquid-like  material  having  excess  free  volume  that 
percolates  through  the  structure  has  been  postulated  earlier  by  Cohen  and  Grest 
[32,33]  in  relation  to  a  glass  transition  and  establishment  of  overall  liquid-like 
behavior  in  bulk.  Our  simulations  are  in  support  of  these  views. 

Our  simulation  has  demonstrated  further  that  glass  transition  takes  place 
when  in  an  undercooled  melt,  the  percolation  condition  established  by  the  strings 
of  liquid-like  material  is  broken.  This  view  is  also  not  entirely  new,  and  while  it 
has  been  set  forth  from  time  to  time  by  others  [see  36],  it  has  not  been  widely 
favored. 

A  related  and  complementary  concept,  that  of  percolation  of  a  property 
called  rigidity,  was  developed  independently  at  the  same  time  by  Phillips  [34]. 
This  was  later  quantified  by  Thorpe  and  formulated  as  a  process  known  as  vec¬ 
tor  percolation  [35].  In  this  concept,  one  speaks  of  random  networks  consisting 
of  floppy  (or  spongy)  and  rigid  regions;  isolated  rings  (polygons)  with  less  than 
6  sides  would  be  rigid,  while  rings  with  7  or  more  sides  would  be  floppy.  When 
rings  become  part  of  a  network,  the  rigid  rings  remain  rigid,  but  a  floppy  ring 
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now  can  become  rigid  if  it  is  next  to  a  rigid  ring,  and  in  this  manner,  rigid¬ 
ity  can  percolate  through  the  network  [35].  Although  rigidity  percolation  was 
developed  for  network  glasses,  and  has  other  applications  to  gellation,  free  vol¬ 
ume  percolation  on  the  other  hand  is  more  pertinent  to  metallic  glasses  and 
the  establishment  of  fluidity.  We  see,  nevertheless,  a  certain  analogy  between 
these  basic  concepts  and  the  percolation  of  structural  dipoles  observed  in  our 
simulation. 

4.2  Topological  Defects 

The  association  of  the  5  and  7  sided  polygons  into  a  structural  dipole  with 
a  definite  increment  of  free  volume  makes  it  a  characteristic  unit  in  the  descrip¬ 
tion  of  the  molten  state,  as  well  as  the  glassy  state  in  two-dimensional  condensed 
matter.  It  is  interesting  that  both  in  the  molten  state  at  the  melting  temper¬ 
ature,  and  in  the  glassy  state  below  Tt,  disorder  in  the  material  is  not  spread 
out  uniformly  over  all  volume  elements.  Rather,  a  considerable  fraction  of  the 
material  at  any  instance  in  time,  maintains  a  reasonable  degree  of  short  range 
order  while  squeezing  the  disorder  into  the  borders  between  the  ordered  regions, 
where  it  is  accommodated  in  the  form  of  the  5-7  sided  structural  dipoles.  Thus, 
this  arrangement  must  have  a  lower  energy  state  than  the  one  where  atomic 
disorder  is  quasi-uniformly  spread  out  over  the  entire  (2-D)  space.  Clearly  also, 
the  5-7  sided  dipole  cannot  be  rearranged  into  a  pair  of  distorted  6  sided  av¬ 
erage  polygons,  since  that  would  have  a  lower  free  volume,  which  apparently 
is  prevented  from  becoming  that  by  the  constraint  imposed  by  the  surrounding 
regions  of  short  range  ordered  material. 
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The  5-7  sided  dipole  of  Voronoi  polygons  combining  a  region  of  excess  di¬ 
latation  and  one  with  a  smaller  compaction  would  appear  to  have  the  features 
of  an  edge  dislocation,  and  could  readily  be  interpreted  as  such.  To  explore 
this  possibility,  we  have  analysed  a  number  of  well-defined  lattice  defects,  such 
as  single  vacancies,  di-vacancies,  tri-vacancies,  and  tetra-vacancies,  edge  dislo¬ 
cations,  and  stacking  faults  that  have  been  found  in  Bragg  bubble  rafts,  which 
have  the  advantage  of  being  readily  available  form  the  literature  [37],  and  being 
fully  equilibrated  through  an  inter-bubble  potential  resembling  closely  certain 
inter-atomic  potentials  of  metals  [18].  Figure  13a-13d  show  respectively  a  va¬ 
cancy,  di-vacancy,  tri-vacancy,  and  tetra-vacancy  in  bubble  rafts,  together  with 
their  Voronoi  polygon  representations.  We  see  that  in  each  case,  these  vacancy 
type  clusters  are  delineated  by  a  pair  of  5-7  sided  dipoles,  in  which  the  7  sided 
Voronoi  polygons  face  each  other,  and  their  interfaces  with  the  5  sided  polygons 
are  separated  by  from  1  to  4  intervening  polygons.  If  the  vacancies  were  aligned 
in  a  line,  such  arrangements  can  indeed  be  represented  by  a  pair  of  edge  dislo¬ 
cations  with  their  tensile  fields  facing  each  other,  terminating  the  rows  of  1,2,3, 
and  4  vacancies.  Figures  14a  and  14b  [37,38]  show  two  edge  dislocations  in  two 
bubble  rafts  with  their  associated  Voronoi  polygon  descriptions.  The  dislocation 
in  Fig.  14a  is  characterized  by  two  adjoining  5  sided  polygons  facing  two  rather 
dilated  6  sided  polygons  across  a  common  straight  border,  while  the  dislocation 
in  Fig.  14b  is  represented  by  a  dipole  of  5  and  7  sided  polygons.  It  is  easy  to 
demonstrate  that  if  the  dislocation  in  Fig.  14a  is  translated  slightly  toward  the 
left  or  the  right,  one  5  sided  polygon  can  be  transformed  into  a  distorted  6  sided 
polygon,  while  the  large  6  sided  polygon  below  the  remaining  5  sided  polygon 
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is  transformed  into  a  7  sided  polygon  to  look  just  like  the  dislocation  in  Fig. 
14b.  Figures  15a  and  15b  show  two  stacking  faults  in  a  bubble  raft  [37].  While 
the  one  in  Fig.  15a  is  planar,  and  is  terminated  by  two  partial  dislocations, 
the  one  in  15b  is  complex  and  distributed  over  a  volume.  Thus,  in  all  of  these 
instances,  in  the  description  of  these  crystal  defects,  edge  dislocations  play  a 
role.  When  the  dipoles  occur  in  closely  spaced  pairs  with  the  7  sided  polygons 
facing  each  other  within  one  or  two  atomic  positions,  the  defects  are  clearly  of 
vacancy  character  in  a  perfect  crystal  environment,  and  an  element  of  large  free 
volume  in  an  amorphous  environment  of  distorted  hexagons.  Such  defect  pairs 
were  found  in  significant  concentrations  in  amorphous  Bragg  bubble  rafts  by 
Argon  and  Shi  [19]  as  sets  of  closely  spaced  Burgers  vector  couples. 

In  the  strings  of  percolating  liquid-like  material  in  the  melt,  subcooled  melt, 
and  the  rapidly  quenched  glass,  the  5-7  sided  dipoles  in  touching  configura¬ 
tions  are  predominantly  ordered  in  sequences  of  5-7-5-7,  etc.  This  indicates  that 
the  material  of  the  strings  acts  in  the  nature  of  a  high  angle  grain  boundary 
with  the  geometrically  necessary  dislocations  being  placed  so  close  as  to  have 
their  cores  penetrate  into  each  other.  Thus,  as  is  well  known  from  the  study  of 
such  high  angle  boundaries,  their  dislocation  nature  should  then  be  completely 
lost,  and  their  topological  contents  must  be  viewed  as  sets  of  polyhedra  (poly¬ 
gons  in  2-D)[39].  Hence,  the  structure  of  the  melt  consists  of  small  islands  of 
quasi-crystalline  order  separated  by  boundaries  of  quasi-periodically  arranged 
polygons  of  5-7  sides,  retaining  their  topological  signatures  as  stacked  edge  dis¬ 
location  cores.  As  we  will  demonstrate  in  (II)  [29],  in  the  melt,  the  liquid-like 
material  (boundary  material)  typically  makes  up  a  fraction  of  about  0.4  of  the 
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total.  Thus,  the  quasi-crystalline  islands  that  are  enclosed  within  the  bound¬ 
aries  are  exceedingly  small  indeed,  and  are  internally  somewhat  disordered  as 
well.  It  is  interesting  to  note  that  the  make-up  of  static  amorphous  Bragg  soap 
bubble  rafts  composed  of  two  dissimilar  size  bubbles  and  studied  by  Argon  and 
co-workers  [8,19]  had  structures  very  similar  to  the  structure  of  the  2-D  mats 
in  the  present  simulations.  In  addition,  the  melting  simulation  carried  out  by 
Fukushima  and  Ookawa  [37]  by  mechanically  vibrating  Bragg  soap  bubble  rafts 
made  up  of  bubbles  of  one  size  have  also  given  melt  structures  that  were  re¬ 
markably  similar  to  results  of  our  simulation,  even  though  in  their  case,  the 
amplitudes  of  vibration  were  largely  normal  to  the  plane  of  the  raft,  rather  than 
being  in  the  plane. 

The  role  of  dislocations  in  melting  and  dislocation  theories  of  melting  are  not 
new  (for  a  discussion  of  some  early  attempts,  see  Nabarro  [40]).  Most  of  such 
developments  were  based  on  the  observation  that  the  energies  of  arrays  of  dis¬ 
locations  increase  less  than  linearly,  while  the  configurational  entropy  increases 
linearly,  resulting  in  a  well-defined  discontinuous  behavior  at  a  definite  tem¬ 
perature  with  a  critical  concentration  of  dislocation.  Among  the  more  modern 
studies  is  the  computer  simulation  of  Koizumi  and  Ninomiya  [41],  who  have  at¬ 
tained  a  glassy  state  by  progressively  increasing  the  density  of  screw  dislocations 
of  all  types  in  a  diamond  cubic  or  f.c.c.  lattice,  followed  by  relaxing  the  positions 
of  the  atoms  on  the  basis  of  an  atomic  pair  potential.  These  authors  have  re¬ 
ported  reaching  RDFs  that  were  quite  similar  to  those  obtained  for  amorphous 
metals,  with  total  densities  of  dislocations  in  the  range  of  1014  -  1018  lines/cm1, 
i.e.,  where  the  dislocation  cores  occupy  nearly  the  entire  volume.  While  these 


authors  do  not  report  the  final  structure  in  more  detail,  they  make  the  obser¬ 
vation  that  core  structures  of  dislocations  remain  intact  in  the  final  amorphous 
material.  Similar  results  were  reported  also  by  Kotze  and  Kuhlmann-Wilsdorf 

[«|. 

The  above  observations  bring  to  fore  long  standing  controversies  on  the  struc¬ 
ture  of  the  amorphous  state  and  more  particularly  the  role  generalized  disloca¬ 
tions  could  play  in  the  inelastic  deformation  of  such  media  (for  a  discussion,  see 
Argon  [43]).  We  will  demonstrate  in  our  associated  simulation  on  plastic  flow 
[31]  (referred  to  earlier  as  (IV))  in  the  two-dimensional  amorphous  media,  that 
while  most  of  the  inelastic  deformation  is  assignable  to  local  shear  transforma¬ 
tions  occurring  inside  the  liquid-like  materials  of  the  boundaries,  a  small  com¬ 
ponent  of  strain  also  comes  from  mobility  of  dislocations  in  the  quasi-crystalline 
domains  encapsulated  by  the  liquid-like  material.  As  can  be  expected,  this  latter 
component  of  deformation  assumes  greater  importance  in  well-relaxed  glasses, 
with  a  high  crystallinity  where  well-formed  and  isolated  dis  locations  are  readily 
recognizable  inside  the  ordered  region,  as  can  be  seen  in  Fig.  16  at  sites  A  and 
B. 

4.3  Generalization  of  Results 


Although  our  observations  of  the  thermo-physical  properties  of  a  two- 
dimensional  idealized  medium  and  its  phase  transitions  qualitatively  reproduce 
nearly  all  features  that  three-dimensional  media  exhibit,  the  results  should  be 
taken  with  some  degree  of  caution.  First,  constraining  the  atoms  to  exist  only 
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in  two-dimensional  mats  prevents  the  local  configurations  from  reaching  lower 
energy  forms  by  moving  in  all  three  directions.  This  invariably  incorporates 
a  higher  free  volume  into  the  structure  than  is  likely  to  be  the  case  in  three- 
dimensions,  and  may  artificially  enrich  the  liquid-like  material.  This  must  be 
one  reason  why  the  melting  points  found  in  our  simulations  were  significantly 
lower  and  the  coefficients  of  expansion  higher  than  expected  in  comparison  with 
three-dimensional  material. 

In  three  dimensions,  we  visualize  the  string-shaped  percolation  of  the  liquid¬ 
like  material  to  take  the  form  of  tortuous  boundaries  pervading  space  and  encap¬ 
sulating  small  quasi-crystalline  domains.  The  projection  into  three  dimensions 
of  closely  spaced  5-7  sided  polygons  representing  liquid-like  material  require  re¬ 
placing  them  with  pairs  of  polyhedra  of  the  type  that  make  up  high  angle  grain 
boundaries  [39].  In  well- relaxed  glassy  materials  with  a  high  degree  of  crys¬ 
tallinity,  dislocation  lines  should  exist,  and  can  be  of  not  only  edge  type,  but 
also  screw  and  of  mixed  nature.  It  is  unlikely,  however,  that  they  play  any  more 
important  role  in  the  description  of  the  structure  or  its  phase  transitions. 
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NOMENCLATURE 


Dimensional  Parameters 


hi 

k 

m0 

Pi 

rH 

To 

Si 

Ua 

Up 

Vi 

x,y 

E0 

H 

V 

S(i)afhig 

T 

nr(*)*w 

c(*)*v 

P 

*(*%» 

Hm) 

n0 

n 


=  Enthalpy  of  atom  i 
=  Boltzmann’s  constant 
=  Mass  per  atom 
=  Pressure  at  atomic  site  i 
=  Distance  between  atoms  i  and  j 
=  Fundamental  length  scale  in  pair  potential 
=  Entropy  of  atom  i 

=  Displacement  parallel  to  cartesian  axis  a 
=  Displacement  parallel  to  cartesian  axis  0 
=  Velocity  of  atom  i 
=  Cartesian  axes  parallel  to  edges 
of  undistorted  simulation  cell 
=  Elastic  constant  defined  at  atomic  site  i 
=  Energy  scaling  factor  -  binding  energy  in  pair  potential 
=  Average  enthalpy  per  atom 
=  Average  potential  energy  per  atom 
=  Elastic  compliance  defined  at  atomic  side  i 
=  Absolute  temperature 

=  Cartesian  axes  in  the  plane  of  simulation  usually  referring 
to  local  atomic  sites 

=  Total  tangential  shear  strain  at  atomic  site  i 
=  Total  tensor  strain  at  atomic  site  i 
=  Density 

=  Tensor  stress  component  at  atomic  site  i 
=  Maximum  shear  stress  (Mohr  circle)  at  atomic  site  i 
=  Energy  potential  between  atom  pairs  i  and  j 
=  Volume  (Voronoi  polygon  area)  per  atom 
in  reference  crystal 

=  Volume  (Voronoi  polygon  area)  per  atom 
in  disordered  solid 
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Derived  Units 


m0 

= 

Unit  of  mass 

= 

Unit  of  length 

= 

Unit  of  energy 

= 

E0/r0,  unit  of  force 

4/ 

= 

yjE0/m0,  unit  of  velocity 

ir 

— 

E0jr\ ,  unit  of  (2-D)  pressure  or  stress 
(generalized  to  E0jr\  in  (3-D)) 

Po 

= 

m/A0,  unit  of  density 

e 

= 

E0/k ,  unit  of  temperature 

= 

r0\/m0/E0,  unit  of  time 

Normalized  Quantities 


p *  =  p/ 7r,  dimensionless  pressure 

v*  =  v/t/,  dimensionless  velocity 

F*  =  $/Ea,  dimensionless  energy 

F*  =  F/<£,  dimensionless  force 

T *  =  T/tf,  dimensionless  temperature 

p*  =  p/ p0,  dimensionless  density 

<7*y  =  <rx y/n,  dimensionless  tensor  stress  element 

t*  =  r/ 7r,  dimensionless  maximum  shear  stress 

fl*  =  n/no,  dimensionless  atomic  volume 


35 


Fig.  1  The  interatomic  pair  potential  between  Cu  atoms  used  in  the  simula¬ 


tion  with  the  two-component  Cu  —  Zr)  material. 

Fig.  2  Temperature  variations  of  dimensionless  volume  fi’,  enthalpy  if*,  and 
potential  energy  V*,  all  on  a  per  atom  basis  for  the  one-component 
model  at  external  pressure  of  p*  =  0.  Heating  and  cooling  results  are 
denoted  by  squares  and  circles  respectively. 

Fig.  3  Same  as  Fig.  2,  except  p*  =  1.0. 

Fig.  4  Same  as  Fig.  3,  except  the  system  is  the  two-component  (Cu  —  Zr) 

model. 

Fig.  5  Temperature  dependent  alterations  in  the  RDF  in  the  two-component 
material  upon  heating  at  p*  =  1.0:  (a)  T*  =  0.05,  (b)  T *  =  0.225, 

(c)  T*  =  0.25  (T  =  Tm),  and  (d)  T'  =  0.3. 

Fig.  6  The  RDF  peak  positions  in  a  perfect  hexagonal  structure. 

Fig.  7  The  temperature  dependent  alterations  in  the  RDF  in  the  one- 
component  material  upon  quenching  at  p*  =  0:  (a)  T*  =  0.35, 

(b)  T*  =  0.25  (T  =  Tm ),  (c)  T*  =  0.15,  and  (d)  T*  =  0.05. 

Fig.  8  Temperature  dependent  alterations  in  the  RDF  in  the  two-component 
material  below  the  melting  point  at  an  external  pressure  of 
p*  =  1.0:  (a)  T*  =  0.24,  (b)  T*  =  0.14,  and  (c)  T*  =  0.06. 

Fig.  9  Temperature  dependent  alterations  in  the  structure  of  the 

one-component  material  upon  heating  at  an  external  pressure 

of  p*  =  1.0:  (a)  T*  =  0.3  (T  =  0.5 Tm),  (b)  T*  =  0.6  (T  =  Tm), 
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incomplete  melting  after  only  10s  time  steps,  and 
(c)  T*  =  0.6  ( T  =  rm),  complete  melting 
after  4000  time  steps. 


Fig.  10  Temperature  dependent  alterations  in  the  structure  of  the 

two-component  material  upon  heating  at  an  external  pressure  of 
p *  =  1.0:  (a)  T*  =  0.1,  the  dark  atoms  are  Zr,  the  light 
ones  Cu,  (b)  T*  =  0.3  (T  =  1.2 Tm)  after  complete 
melting  the  5  and  7  sided  polygons  associate  into  percolating 
strings  that  traverse  across  the  cell. 

Fig.  11  Temperature  dependent  alterations  in  the  structure  of  the 

two-component  material  upon  quenching  at  an  external  pressure 

of  p*  =  1.0:  (a)  T*  =  0.3  (T  =  Tm),  (b)  T*  =  0.2 

(T,  <  T  <  Tm ),  (c)  r*  =  0.14  (T  <  r,),  and  (d)  T*  =  0.06. 

Fig.  12  Dependence  of  the  average  volume  of  the  Voronoi  polygons 

on  the  number  of  sides  of  the  polygons,  in  the  two-component 
material  at  T*  =  0.1. 

Fig.  13  Distortions  produced  in  the  Voronoi  polygon  field  of  a  perfect 
hexagonal  buble  raft  produced  by  specific  vacancy  clusters: 

(a)  vacancy,  (b)  di- vacancy,  (c)  tri-vacancy,  and 
(d)  tetra-vacancy. 

Fig.  14  Distortions  produced  in  the  Voronoi  polygon  field  of  a  perfect 
hexagonal  buble  raft  produced  by  edge  dislocations:  (a)  edge 
dislocation  in  an  intermediate  position  between  two  low  energy 
configurations,  and  (b)  edge  dislocation  in  the  low  energy 
configuration. 
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Fig.  15  (a)  A  stacking  fault  with  two  terminating  partials  in  a  perfect 

hexagonal  bubble  raft,  and  (b)  a  complex  two-dimensional  fault. 

Fig.  16  Some  isolated  edge  dislocations  (5-7  sided  dipoles)  in  a  well 
relaxed  two-component  glass. 
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Temperature  dependent  changes  in  the  volume  per  atom  ft*,  en 
thalpy  per  atom  H*,  and  potential  energy  per  atom  V*,  all  i 
normalized  units  for  the  one-component  material  under  an  ex 
temal  pressure  of  p  =  1.0. 
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Temperature  dependent  alterations  in^the  RDF  in  the  two- 
component  material  ugon  heating  at  p  =  1.0:  (a)  T *  0.05, 

(b)  T*=  0.225,  (c)  r-  0.25  (T  =  Tm),  and  (d)  T  =  0.3. 
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Temperature  dependent  alterations  in  the  RDF  in  the  two- 
component  material  below  the  melting  point  at  an  external 
pressure  of  p  =  1.0:  (a)  T  =  0.24,  (b)  T  -  0.14,  and 

(c)  ?•  0.06. 


Temperature  dependent  alterations  in  the  structure  of  the  one- 
component  material  upon  heating  at  an  external  pressure  of 
p*=  1.0:  (a)  r»  0.3  (T  =  0.5  T  ),  (b)  T  =  0.5, 

(c)  T  =  0.6  (T  %Tm)>  incomplete  melting  after  only  103  time 
steps,  and  (d)  T  -0.6  (T  s  T  )  complete  melting  after  4000 
time  steps.  m 
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Fig. 9c 


Temperature  dependent  alterations  in  the  structure  of  the  two- 
component  material  upon  heating  at  an  external  pressure  of 
p*=  1.0:  (a)  T  = =  0.1, ^the  dark  atoms  are  Cu,  the  light  ones 

Zr,  (b)  T  =  0.3  (T  -  T  )  after  complete  melting 

the  5  and  7  sided  polygons  associate  into  percolating  strings 
that  traverse  across  the  cell. 
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Fig.  11  -  Temperature  dependent  alterations  in  the  structure  of  the  two- 

component  material  upon  quenching  at  an  external  pressure  of 
P*=  1.0:  (a)  T*=  0.3  (T  =  Tm) ,  (b)  f  =  0.2  (T  <  T  <  T  J , 

(c)  f -  0.14  (T  <  T  ) ,  and  (d)  0.06.  9 
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AVE.  AREA  OF  V.P. 


SIDE  NUMBER  OF  V.P. 


Dependence  of  the  average  volume  of  the  Voronoi  polygons 
the  number  of  sides  of  the  polygons,  in  the  two-component 
material  at  T  -  0.1. 
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Fig.  13  -  Distortions  produced  in  the  Voronoi  polygon  field  of  a  perfect 

hexagonal  bubble  raft  produced  by  specific  vacancy  clusters: 
(a)  vacancy,  (b)  di-vacancy,  (c)  tri-vacancy,  and  (d)  quadri- 
vacancy. 
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Distortions  produced  in  the  Voronoi  polygon  field  of  a  perfect 
hexagonal  bubble  raft  produced  by  edge  dislocations:  (a)  edge 
dislocation  in  an  intermediate  position  between  two  low  energy 
con  iguraticrs,  ana  vb)  edge  dislocation  in  the  low  energy  con¬ 
figuration. 


(a)  A  stacking  fault  with  two  terminating  partials  in  a  perfect 
hexagonal  bubble  raft ,  and  (b)  a  complex  two-dimensional  fault. 
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Fig.  7 5b 


Some  isolated  edge  dislocations  (5-7  sided  dipoles)  in  a  well 
relaxed  two-component  glass. 


